BigDFT.Interop.ASEInterop module

This module contains some wrappers for using ASE to perform various operations on BigDFT molecules.

https://wiki.fysik.dtu.dk/ase/

class BigASECalculator(parameters, calc, restart=None, ignore_bad_restart_file=False, directory='.', label=None, atoms=None, **kwargs)[source]

Bases: ase.calculators.calculator.Calculator

This calculator can be used by ASE to run BigDFT calculations.

Parameters
calculation_required(atoms, quantities)[source]

https://wiki.fysik.dtu.dk/ase/_modules/ase/calculators/interface.html

get_forces(atoms)[source]

https://wiki.fysik.dtu.dk/ase/_modules/ase/calculators/interface.html

get_potential_energy(atoms=None, force_consistent=False)[source]

https://wiki.fysik.dtu.dk/ase/_modules/ase/calculators/interface.html

get_stress(atoms=None)[source]
ase_potential_energy(sys, ase_calculator)[source]

Given a ASE calculator, calculates the potential energy of a BigDFT System.

Parameters
  • sys (BigDFT.Systems.System) – a BigDFT system

  • ase_calculator (ase.calculators.calculator.Calculator) – ASE calculator

Returns

the potential energy, in Hartree

Return type

float

ase_to_bigdft(mol)[source]

Given an ASE collection of atoms, this transforms that collection into a BigDFT fragment.

Parameters

mol (ase.Atoms) – a collection of atoms used by ASE.

Returns

a BigDFT fragment.

Return type

(BigDFT.Fragments.Fragment)

bigdft_to_ase(sys)[source]

Given a BigDFT system, this transforms that system into an ASE collection of atoms.

Parameters

frag (BigDFT.Systems.System) – the system to convert.

Returns

a collection of atoms for ASE.

Return type

(ase.Atoms)

BigDFT.Interop.ASEInterop example

Below we show an example of using the ASE interoperability functions.

def _example():
    """Example of using ASE interoperability"""
    from BigDFT.IO import XYZReader
    from BigDFT.Inputfiles import Inputfile
    from BigDFT.Calculators import SystemCalculator
    from BigDFT.Fragments import Fragment
    from BigDFT.Systems import System
    from ase.calculators.lj import LennardJones
    from BigDFT.UnitCells import UnitCell
    from ase.units import Hartree
    from ase.optimize import BFGS
    from os.path import join
    from copy import deepcopy

    # Create a system.
    sys = System()
    reader = XYZReader("Ar")
    sys["FRA:1"] = Fragment(xyzfile=reader)
    sys["FRA:2"] = deepcopy(sys["FRA:1"])
    sys["FRA:2"].translate([-2, 0, 0])

    # Skip straight to the potential energy.
    print(ase_potential_energy(sys, LennardJones()))

    # Advanced used.
    asys = bigdft_to_ase(sys)
    asys.set_calculator(LennardJones())
    dyn = BFGS(asys)
    dyn.run(fmax=0.05)
    print(asys.get_potential_energy() / Hartree)

    # Unit cells
    sys.cell = UnitCell([5, 5, 5], units="bohr")
    print(ase_potential_energy(sys, LennardJones()))
    sys.cell = UnitCell([5, float("inf"), 5], units="bohr")
    print(ase_potential_energy(sys, LennardJones()))
    sys.cell = UnitCell([float("inf"), float("inf"), 5], units="bohr")
    print(ase_potential_energy(sys, LennardJones()))

    # We can also use BigDFT with ase.
    inp = Inputfile()
    inp.set_xc("PBE")
    inp.set_hgrid(0.4)
    code = SystemCalculator(verbose=False)
    sys.cell = UnitCell()
    print(ase_potential_energy(sys, BigASECalculator(inp, code,
                                                     directory="work",
                                                     label="ase-free")))
    sys.cell = UnitCell([5, 5, 5], units="bohr")
    print(ase_potential_energy(sys, BigASECalculator(inp, code,
                                                     directory="work",
                                                     label="ase-periodic")))
    sys.cell = UnitCell([5, float("inf"), 5], units="bohr")
    print(ase_potential_energy(sys, BigASECalculator(inp, code,
                                                     directory="work",
                                                     label="ase-surface")))
    sys.cell = UnitCell([float("inf"), float("inf"), 5], units="bohr")
    print(ase_potential_energy(sys, BigASECalculator(inp, code,
                                                     directory="work",
                                                     label="ase-wire")))